
use $datapath/clean/main_data.dta, clear 

************************************************
* Create variable labels: 

	label var Qp "phys"
	label var Np "phys holdings"
	label var Qe "ebook"
	label var Ne "ebook holdings"

			
**********************************************************************
** create jump variables


************************************
** ebook jumps 
	
	gsort lno year 
	
** levels: 	
	gen ex = Ne-Ne[_n-1] if lno==lno[_n-1] 
	su ex, de
		local etd = r(p90)

** jump cutoff: 
gsort lno year
	gen ee = (Ne-Ne[_n-1] >`etd') & lno==lno[_n-1] & (Ne+Ne[_n-1]~=.) & (Ne[_n+1]-Ne >=0)
	
******************


******************
** create IVs

	gsort lno year
		by lno: gen ive =sum(ee) 
		label var ive "ebook jump"

** First jump year:
	gen xe = year if ive==1  
	egen Ymin = min(xe), by(lno)

** total number of jumps:  
	egen total_jumps = sum(ee), by(lno)	

** Set up observations we want to include 
	gen ifstatement = (ive <2 & year<=Ymin+1 & year >= Ymin-2) | total_jumps==0
	
**********************************************************************


		
**********************************************************************
** run regressions & create table:
	
	xtset lno year
	eststo clear
		
** OLS effect on Qp & Qe
	eststo: reg visits Np Ne i.year, vce(cluster lno)

** FE effect on Qp & Qe
	eststo: areg visits Np Ne i.year, absorb(lno) vce(cluster lno)
		estimates save output/causal_visits, replace

	
** For text: reduced form
	areg visits ive i.year if ifstatement==1, absorb(lno) vce(cluster lno)

		
** IV effect on Qp & Qe
	eststo: xtivreg visits Np ( Ne = ive) i.year if ifstatement==1, fe vce(cluster lno)	

	
	
** create table 

	local pattern prefix(/multicolumn{@span}{c}{) suffix(}) span erepeat(/cmidrule(lr){@span})

	esttab,  noomitted mtitles("OLS" "FE" "IV") drop(*.year _cons) se replace  label nonotes star(* 0.10 ** 0.05 *** 0.01) scalars("r2_a $/overline{R^2}$") nocons

	esttab using "$tablepath/table4_visits_causal.tex",  noomitted mtitles("OLS" "FE" "IV") drop(*.year _cons) se replace  label nonotes star(* 0.10 ** 0.05 *** 0.01) scalars("r2_a $/overline{R^2}$") nocons booktabs 
**********************************************************************






**********************************************************************
** numbers in text 


************************
** decline in visits due to ebooks existing

areg visits Np Ne i.year, absorb(lno) vce(robust)
	gen bfe=_b[Ne]
	su Ne  if year==2018
	local m=r(mean)*r(N) 

	su visits if year==2018
	local v=r(mean)*r(N) 

	gen decline_fe = `m'*bfe/`v'

************************
** decline in visits between 2013 and 2018
	
	areg visits i.year, absorb(lno)
		local dv = _b[2018.year]
		
	areg Ne i.year, absorb(lno)
		local dn = _b[2018.year]
		
	gen decline_1318 = `dn'*bfe / `dv'
	su decline*
	
	
************************
** Decline in visits

	areg visits i.year, absorb(lno)
		di _b[2018.year] / _b[_cons]







